Microscopic Calculation of Flow Stress in Cu-Mg Metallic Glass 
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We have carried out shear-deformation simulations on amorphous Mg-Cu systems at zero temper- 
ature and pressure, containing 2048-131072 atoms. At the largest size a smooth stress-strain curve 
is obtained with a well-defined flow stress. In the smallest system there are severe discontinuities in 
the stress-strain curve caused by localized plastic events. We show that the events can be charac- 
terized by a slip volume and a critical stress and we determine the distribution of these quantities 
from the ensemble of all events occurring in the small system. The distribution of critical stresses 
at which the enthalpy barriers for the individual events vanish is spread between 200 MPa and 500 
MPa with a mean of 316 MPa, close to the flow stress observed in the largest system. 



The mechanical properties of bulk metallic glasses [l|,|2| 
(BMGs) are the subject of intense research. A host of 
applications is envisaged if only reasonable macroscopic 
plasticity could be achieved, rather than the intense lo- 
calization into shear bands which typically occurs. De- 
tailed knowledge of plastic deformation mechanisms in 
glasses, and their connection to macroscopic flow prop- 
erties, however, remains elusive: While in crystals the 
dislocation provides a well defined starting point for es- 
timates of flow stress, in glasses there is no such easily 
characterizable defect. 

Various recent theories Q, 0] take as a starting point 
the existence of a collection of "relaxation centers" [3 or 
"shear-transformation zones" (STZs)|j,|5| which operate 
as localized centers of deformation. Based on a set of 
assumptions about the operation of the relaxation centers 
the total rate of plastic deformation is obtained from the 
collective behavior of the local centers. 

Indeed, several simulations of deformation in amor- 
phous metals la, ly, LE 1^ 1^ have established that the 
plastic behavior involves localized events with up to 150 
atoms 6] . The transforming regions have been exten- 
sively studied from a p otential energy landscape point 
of view by Lacks [a IllJ, luj, and the response to shear 
and normal stresses have been investigated for Lennard- 
Jones systems [l^- In the latter work the connection be- 
tween the properties of a single model STZ and the yield 
stress was discussed through the application of the Mohr- 
Coulomb yield criterion. Systematic studies of the depen- 
dence of macroscopic properties on the full ensemble of 
local deformation events, however, are still needed. 

In this Letter, we describe zero temperature simula- 
tions of plastic deformation of a Cu-Mg glass and demon- 
strate a direct connection between the statistical proper- 
ties of localized deformation events and the flow stress 
at mesoscopic length scales (10-100 nm). Our main 
results are (i) while for small systems the stress-strain 
curve shows sharp drops signifying individual deforma- 
tion events, when the system size is over 10^ atoms, the 



stress strain curve becomes quite smooth, with a clear 
flow stress of 320-330MPa (FigP); (ii) analysis of the in- 
dividual events in the smallest system can be quantified 
in terms of two quantities, a critical stress Uc and a quan- 
tity we call the "slip volume" V^np; (iii) the mesoscopic 
flow stress is the mean of the distribution of Uc, and the 
spatial density of transforming regions per unit strain in 
a large system is the inverse of the mean of Kiiip • 

The simulated material is Mgo.85Cuo.15, which is 
the optimal glass-forming composition for the Mg-Cu 
system^lS^. This system is interesting because the 
addition of a small amount of Y makes it a BMG. 
The interatomic potential is the effective medium the- 
ory (EMT)[lJ|, fitted to properties of the pure ele- 
ments and intermetallic compounds obtained from ex- 
periment and density functional theory calculations. De- 
tails of the potential and of the method for creating the 
zero-temperature glassy configurations may be found in 
Ref.USt for the 16384- and 131072-atom systems we have 
used the cooling rate of 0.72 K/ps. The main simu- 
lations involve continuous deformation of systems con- 
taining 2048, 16384 and 131072 atoms. The nominally 
zero-temperature configurations resulting from the cool- 
ing runs were further minimized before deformation runs, 
with respect to both atomic positions and the vectors de- 
scribing the periodic supercell. The box vectors are in 
fact controlled by six strain degrees of freedom, which 
also play a role in the deformation simulations. We re- 
strict to pure shear. The deformation simulations are 
strain-controlled: the relevant component of strain is in- 
cremented in steps of 0.0005. After each step the remain- 
ing degrees of freedom — the atomic positions and other 
strain components — are relaxed to minimize the energy. 

Fig- n shows stress-strain curves for shear deformation 
simulations for three different system sizes. In all cases 
there is initially a smooth, almost linear increase, mark- 
ing the elastic regime. In the smallest system, this is 
first interrupted by small kinks in the stress-strain curve, 
and then by sharp drops in stress. The remainder of the 
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FIG. 1: Upper panel, stress-strain curves for different-sized 
systems undergoing pure shear at zero temperature, shifted 
for clarity. For the largest system the peak and flow stresses 
are 450MPa and 320MPa respectively, or about /i/16 and 
^/23 where /i=7.3GPa is the shear modulus. Lower panel, 
stress-strain curve of a sub-volume of the 131072-atom sys- 
tem, selected because visualization (via deviation from afflne 
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FIG. 2: Distribution of CTc and VsUp values for all events (un- 
shaded) and events taking place above 10% strain (shaded). 
The means are 331(316) MPa and 305(317) A^, where the 
flgures in parentheses indicate the means of the shaded dis- 
tributions. 



curve is characterized by intervals of almost linear in- 
crease punctuated by further sharp drops. The highest 
stress attained is prior to the first significant drops. The 
size of fluctuations (drops) progressively diminishes for 
larger system sizes, so that for the 131072-atom system 
we observe a quite smooth curve with a well-defined flow 
stress following an initial peak. 

A natural interpretation of the stress drops in the 
small system is that they represent discrete plastic flow 
"events" , corresponding to the postulated STZs of the- 
ory or the localized events found in previous simulations. 



The fact that the slope of the stress curve between the 
drops is always the same indicates that here only elastic 
deformation is taking place and that the plastic defor- 
mation is entirely accounted for by processes associated 
with the drops. In the following we shall analyze the 
individual events, as identified in the small systems and 
use this information to address the properties of the large 
system. To confirm that the deformation of the large sys- 
tem also happens through localized events, we have used 
the technique from Ref. |a to highlight atoms with a large 
deviation from affine deformation |2l| . and observing by 
direct visualization that the atoms so highlighted tend 
to form clusters, have generated statistics of these clus- 
ters. A cluster is a group of (at least three) so-highlighted 
atoms connected by nearest neighbor bonds. We count 
5-7 clusters per nm'^ per unit strain in the "steady-state" 
regime between 10% and 20% strain. These clusters typi- 
cally contain 3-20 atoms, although extreme cases involv- 
ing up to 125 atoms also occur. In the small system an 
event marked by a stress drop typically involves one or 
two such clusters, thus we can use this system to study 
single events in detail. 

The events are transitions involving internal rearrange- 
ments of atoms. We have identified two characteristic 
quantities associated with events. The first of these, 
termed "slip volume" , is geometrical in nature and rep- 
resents the amount of plastic strain associated with the 
event. The plastic strain is defined in terms of changes in 
the shape of the periodic simulation cell. However, if we 
take the event to be localized in the cell, the plastic strain 
induced at the boundaries depends not just on the geom- 
etry of the rearranging atoms but also on the volume of 
the cell. To see this, consider an idealized slip event as a 
planar area Ah cut within the material, and the resulting 
free surfaces shifted relatively by an amount 6 _L n, as in 
dislocation loop nucleation. Then the plastic shear strain 
felt by the boundary of the system is epi ^ Kiiip/1^, where 
V^iip — bA, and V is the system volume. Alternatively, 
knowing e^/ we can multiply by V to obtain V^up. This 
slip volume is a tensorial quantity, but in this work we 
are only concerned with the component corresponding to 
the applied shear strain. We cannot decompose it into 
b and An, since the events are in general geometrically 
more complex than planar slip. 

As we describe below we have determined the distribu- 
tion of 14iip, shown in the right panel of Fig. |2 Rather 
than being peaked at a finite value, the most likely value 
tends to zero. The distribution is in fact roughly expo- 
nential, although with some extra weight at large values: 
The mean is t^up ^ 305 A'^; if an exponential is fitted to 
the initial data, a slightly smaller characteristic T^iip of 
240 A^ is obtained. In a large system, the total strain 
can be written etot = NeVsUp/V, where N^ is the total 
number of events. This implies that V^- = 3.2nm^'^ is 
the number density of events per unit strain and volume. 
This is roughly a factor of two smaller than the num- 



ber obtained by counting clusters; this is partly due to 
the fact that in the small system energetically isolated 
events may have more than one spatially separate clus- 
ter, and partly due to a possible enhancement of spatially 
separate events "cascading" together due to interactions 
with periodic images in a small system. Cascading has 
been studied by Maloney 16]. If we use the estimation of 
the mean from the exponential fit (which ignores excess 
large events) instead of the actual mean, we find an event 
density 4.2nm~'^, which is closer to that determined ge- 
ometrically. 

The second quantity is a critical stress, cTc, the value 
of the shear stress at which the event happens sponta- 
neously at T=0. At lower stresses, there exists a bar- 
rier which might be crossed due to thermal fluctuations 
at finite temperature, but at T=0 prevents the event 
from taking place until the stress rises high enough. To 
properly define Uc we adopt a stress-controlled formal- 
ism, where the six strains all become degrees of free- 
dom, and a stress term is added to the potential energy, 
with stress as a tunable parameter. There then exists a 
(3N-|-6)-dimensional potential energy (or enthalpy) land- 
scape, which is effectively tilted by increasing the stress. 

For each event apparent in the stress-strain curve, we 
take configurations from the simulation before and after 
the stress drop. These are close to minima of the en- 
thalpy landscape. By minimization under a chosen stress 
we obtain locally stable configurations which we take as 
the "initial" and "final" states for the given event and 
the given stress. These are used to define V^up. The 
enthalpy barrier between the events is determined using 
the Nudged Elastic Band method[i3,ElIi3- These are 
computed for a range of stresses sufficient to determine 
tJc with accuracy; they are stress dependent, but the de- 
pendence is small in the case of T4iip, and was averaged 
over. We have taken pairs of configurations and calcu- 
lated the above quantities in this way for every peak on 
the stress-strain curve from a shear deformation simula- 
tion of a 2048-atom system up to 30% strain. 

Enthalpy profiles for a particular event are shown in 
Fig. 121 In the regime where the calculations have been 
done, i.e., fairly close to the critical stress at which the 
barrier vanishes (necessary for events to occur at T=0), 
the barrier is very small (1-10 meV) compared to the 
overall enthalpy change (~1 eV). Also shown in the fig- 
ure are the atoms most involved (determined by deviation 
from afiine deformation). An animation of the process 
shows that the nature and order of the individual mo- 
tions is as roughly indicated by the arrows. This event is 
relatively simple; some involve several tens of atoms and 
complex patterns of motion, although much of this can 
be decomposed into such small pieces involving snake- like 
motion and rotations of groups of three or four atoms. 

The inset of Fig. El shows the barrier's stress depen- 
dence, whose form is a steady decrease with increasing 
stress, flattening out somewhat at the critical stress <Tc 
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FIG. 3: Enthalpy along several minimum energy paths for 
an event, for different applied stresses as indicated; the bar- 
rier is just visible as a slight bump at the beginning of the 
path. In the bottom are shown the participating atoms, col- 
ored according to deviation from afhne deformation computed 
between initial and final states, with arrows giving an indi- 
cation of the types of motion and the numbers identifying 
the order of motions and the corresponding features in the 
enthalpy curve. Inset: stress-dependence of the barrier. 



where the barrier vanishes. In fact, the barrier must van- 
ish with a 3/2 power law sufficiently close to ac-, due to 
the merging of a saddle-point and a local minimum of 
the enthalpy (as in a saddle-node bifurcation). The res- 
olution of our data is not enough to identify this; simply 
using a linear or quadratic fit to the data near Cc is suf- 
ficient to identify CTc with an accuracy of ^ 2 MPa. The 
vanishing of the barrier (and the minimum) has been 
studied in detail by Malandro and Lacks,8,, 10]. At lower 
stresses, there is extra structure in the barrier height, 
such as abrupt changes of slope or even local maxima. 

In some cases an intermediate minimum was found 
along the minimum energy path, and in these cases sepa- 
rate calculations were made for each of the thus-identified 
"sub-events" . By combining the results from all bar- 
rier determinations we can plot the distribution g{(Tc) 
of critical stresses, shown in Fig. (21 left. To improve 
the statistics, we have also included events obtained by 
shear-deformation in the x — z and y — z planes, yielding 
a total of 262 events. There is a broad peak, with a mean 
of 331 MPa and a standard deviation of about 70 MPa. 
If we count only events taking place after 10% deforma- 
tion (the shaded distribution in Fig. O, the mean is a 
little lower, 316 MPa. This latter value is close to the 
flow stress observed in the large system. 

Why is this? To obtain a connection between the de- 
formation behavior of the small system and that of the 
large, we have computed the stress averaged over a sub- 
set of the 131072-atom system whose volume is that of 
the 2048-atom system. This stress-strain curve is plot- 



ted in the lower part of Fig. ^ The striking feature of 
this curve is that it looks closer to that of the actual 
2048-atoni system than it does to that of the 131072- 
atom system as a whole, although the stress drops are 
not quite as large, nor as sharp. The particular subset 
was chosen to surround a cluster which was active at 
strain 0.064, and indeed a significant drop in the stress 
can be seen at that strain. The large system thus in 
a sense behaves as a collection of weakly coupled small 
systems, each undergoing relatively large stress fluctua- 
tions, but whose average is quite smooth. The fact that 
the stress drops are gentler than in the true 2048-atom 
system is presumably due the smaller constraints on this 
region provided by the surrounding material, compared 
to those of periodic boundary conditions; stress can re- 
lax into neighboring material during the relaxation to 
mechanical equilibrium. 

Simple considerations provide crude estimates for the 
mesoscopic flow stress. In the manner of averaging elas- 
tic constants in Dolvcrvstalsj20j. we can assume that ei- 
ther the stress or the strain is uniform over subsystems. 
In the first case we imagine imposing a fixed stress on 
all subsystems and letting them respond independently. 
Then under relaxation, every subsystem will flow until 
it reaches the first critical stress that is higher than the 
imposed stress. No further deformation can take place — 
unless the imposed stress is higher than the maximum 
critical stress. Thus the flow stress is the maximum of 
the (Tc-distribution. This is at least an upper bound: The 
material clearly cannot sustain a stress greater than the 
maximum CTc, 450 MPa for events taking place above 10% 
strain. 

Alternatively, imposing a uniform strain on each sub- 
system, each undergoes deformation just like the sin- 
gle simulated small system. Assuming the individual 
stress-strain curves have no fixed phase relation, aver- 
aging across them at a given strain is equivalent to av- 
eraging over the strain history of a single subsystem. 
This average is straightforward to compute if the assump- 
tion of a fixed T4iip is made, yielding CTc — At7/2, where 
Act = ^iiVsiip/Vsub is the stress drop associated with in- 
dividual events in a subsystem. This value is 265 MPa 
with our choice of small system {Vsub = 4.3 x lO^A^), 
about 10% less than the observed flow stress. 

This estimate is necessarily rough, in particular be- 
cause it explicitly involves an apparently arbitrary sub- 
system size. However, the size of our small system is 
not so arbitrary. We have previously noted that it cor- 
responds more or less to the size at which events be- 
come discrete. Furthermore, analysis of the distribution 
of stress (not shown here) averaged over various-sized 
subsystems suggests that this size is about the smallest 
at which random fluctuations coming from the atomic 



stresses begin to cancel out enough to make the aver- 
aged stress meaningful — a quantity which actually rep- 
resents force per unit area exerted by the material on 
itself. Apart from these considerations, the directly ob- 
served correspondence between the flow stress and the 
mean Cc itself supports the overall picture of the statistics 
of small subsystems determining the mesoscopic plastic 
behavior. 
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